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TECHNICAL MEMORANDUM X- 64833 


ASYMPTOTIC SOLUTION TO THE TANGENTIAL LOW THRUST 
ENERGY INCREASE TRAJECTORY 

1. INTRODUCTION 


A spacecraft, equipped with a low-thrust constant acceleration propulsion device, 
ascending from an initial circular orbit around a spherical central body describes a spiral 
orbit. To determine the radius vector and velocity of the spacecraft at any time, the 
equations of motion need to be integrated numerically, even in the case of a tangential 
thrust steering program which is considered in this study. It is well-known [ 1 ] that for 
energy increase trajectories, the tangential thrust steering is very close to the optimal 
steering program, especially in the inner, multirevolution part of the spiral trajectory, 
which is characterized by the acting of a small force on a space vehicle in the presence of a 
strong gravitational field. With the absence of any other perturbing force, nearly Keplerian 
orbit conditions exist locally. The deviations between the optimal and the tangential thrust 
trajectories arise in the outer part on which the vehicle moves from near circular conditions 
to parabolic velocity and the gravitational force drops to the same order of magnitude as 
the thrust. 

To gain some experience in performing the analysis of a very accurate analytical 
approximate solution of the optimal trajectory, the tangential thrust trajectory offers a 
very good simplified problem, since it possesses all the main characteristics of the optimal 
solution. 

In the literature, a number of analytical approximations of the optimal trajectory 
are presented, and many are available that consider the tangential thrust trajectory [2,3] . 
The method used herein to derive a solution for the tangential thrust trajectory follows 
those of Reference 2. The basic method is described by Cole and Kevorkian [4] based on 
previous work of Linstedt and Poincare [ 5,6] . However, formulation of the equation of 
motion using regularized variables, as presented in Reference 1 , allows one to derive an 
extremely accurate second order solution using a very simple analysis. 

As will be shown in Section IV, in comparison to numerically calculated tangential 
thrust trajectories, the analytical solution derived is extremely accurate. The accuracy 
decreases with time, which indicates that in addition to the basic characteristic limitations 
previously discussed, accuracy requirements will determine the limits of the applicability 
of the solution presented. 

For e < 0.001 , the solution for all of the state variables gives at least three-digit 
accuracy, up to an energy level close to escape. The solution fails to describe the trajectory 
near escape. 



Sections II and III contain the derivation of the solution which is summarized in 
Section III.E. Section IV shows an alternate time approximation. The solution is discussed 
in Section V. 


II. MATHEMATICAL FORMULATION 


The ascent of a spacecraft is considered from an initial circular orbit around a 
spherical planet. The vehicle is described by a point mass and is equipped with a low-thrust 
constant acceleration engine. The engine produces a small continuous thrusting force. 

No other perturbations are present. The trajectory is considered as planar. 

The following figure shows the definitions of symbols used for the mathematical 
description of the spacecraft motion. 



GRAVITATIONAL FORCE CENTER 


The equations of motion are formulated as follows: 


r' = u 

i p = v/r - (1) 

u' = m + 2hr + er 2 - 

|v| 
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V 


V 


( 1 ) 

(Concluded) 

with the absolute value of the vehicle velocity 

|V| = yfu 2 + V 2 , (2) 

and the thrust to mass ratio 

e = F/m . (3) 

The acceleration is considered as constant, 

e = constant . (4) 

Using the initial orbit as a scaling system, the initial conditions are 

r(s = 0) = r(0) = 1 
V>(s = 0) = *>(0) = 0 
u(s = 0) = u(0) = 0 

(5) 

v(s = 0) = v(0) = 1 
h(s = 0) = h(0) = - 0.5 
t(s = 0) = t(0) = 0 

Equation system (1) together with definition equations (2) and (3) and the boundary 
conditions (5) describe the tangential thrust spiral orbit of a low-thrust space vehicle. In 
this particular formulation, the energy h is used as an additional (or redundant) variable. 



h' = e |v| 
t' = r 
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The independent variable s is defined by the last equation of system (1 ). (For detailed 
information, see Reference 1 .) 

Although the system (1 ) cannot be solved analytically, the solution can be approxi- 
mated by use of a two-variable expansion procedure. The variables will be evaluated by an 
asymptotic expansion in powers of the small parameter e, 

x = x (0) + ex* 1 * + e 2 x (2) + e 3 x^ + 0(e 4 ) , (6) 


with x in equation (6) indicating any variable. 

In addition, the variables are considered as functions of two new independent 
variables defined by 


r 


e 


es 


s 

r = f co(r € )do 
0 


(7) 


The derivatives of equations (7) with respect to the independent variables are 


r e = e 


T e ^ = 0 


t' = cj(tJ 


r (n) _ e n 


d”- 1 ^ 


dr, 


n-1 


( 8 ) 


with n > 2. The derivations of a variable with respect to the new independent variables 
follow therefrom: 
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( 9 ) 


x' = x,T e ' + x 2 t' = cx i + co(r e )x 2 

x" = x u (T e ') 2 + x^Ct') 2 + x 2 r'' = 6 2 x j[ + gj 2 x 22 + e 2 — x 2 

37 V 


The subscripts 1 and 2 denote the partial derivatives defined as follows: 



x 


2 


9x 

9r 



= Oh 

~ dr 2 


_ 9 2 x 

9r e 9r 


( 10 ) 


To derive an asymptotic solution of equations (1), the independent variable of the 
system is transformed to the two independent variables defined by equations (7). After 
expanding the right side of equations (1 ) using power series of the variables as given in 
equation (6), collecting terms of the same power of e will lead to sets of partial differential 
equations whose solution is given in Section III. 


III. ASYMPTOTIC EXPANSION SOLUTION 


A. Base or Zero Order Solution 


Performing as stated in Section II will lead to the following zero order differential 
equations: 


5 



cor 2 (0) = 

u (o) 


uv 2 fe > = 

v «»/r 

(0) 

cou 2 to) = 

1 + 

2h (0) r (0) 

cov 2 (°) = 

0 


coh 2 (0) = 

0 


cot 2 (0 > = 

r (o) 



with boundary conditions 


B 

O 

o 

II 

1 

^ (0) (0,0) = 

0 

u (o) (0,0) = 

0 

v (o) (0,0) = 

1 

h (o) (0,0) = 

-0.5 

t (o) (0,0) = 

0 


(ID 


( 12 ) 


Instead of using the complete solution of systems (1 1) and (12) describing a 
generalized Keplerian solution of the equations of motion without thrust, a partial solution 
is used as a base solution as follows: 

r (o) = R«»( Te ) 

= or-' V (0) R (0)1 t + 4> (0 \t € ) (13) 

u (0) = 0 
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v (0 > = v (0) (r e ) 
h (0) = H (0) (r e ) 
t (0) = R (0) co' 1 r + T< 0 >(Tg) 

Inserting this solution into equations (11) will give 

H (0) = - [ 2 R (0) ] _ * . (14) 


(13) 

(Concluded) 


The capital letters denote unknown functions of the slow variable r e only; i.e., constants 

with respect to an integration over the fast variable r. They are determinable by removing 
secular r-terms due to the asymptotic expansion from the first order solution and by satis- 
fying the boundary conditions (12). 


B. First Order Solution 

The first order differential equations derived are 


cor 2 (1) = - R, (0) + u (1) 

uxp 2 (*> = . ^j«» + V (1) /R (0) - V (0) r (1) /R (0)2 
wu 2 (1) = 2 [R (0) h (1) + H (0) r (1) ] 
wv 2 (1 > = - V, (0) + R (0)2 
coh 2 (1) = - Hj (0) + V (0) 


(15) 


wt 2 (1 > = - t, (0 > + r (1) 


with boundary conditions 


(16) 


r (I) (0,0) = 0 

^ ( 1 ) ( 0 , 0 ) = 0 
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u (1) (0,0) = 0 


v (1) (0,0) = 0 
h (1) (0,0) = 0 


t (1) (0,0) = 0 


( 16 ) 

(Concluded) 


Because the right side of the velocity component v and energy equations are independent 
of r, it follows that 

V (D = [. V, (0) + R (0)2 ] t + V (l) (r e ) (17) 

and 

h (1) = [- Hj (0) + V (0) ] t + H (1) (r e ) . (18) 

To determine H (0) , R f0) , and V (0) , the secular terms in equations ( 1 7) and (18) will 
be removed: 

- Vi (0) + R (0)2 =0 (19) 

and 

- H, (0) + V (0) = 0 . (20) 

Using equation (14) we find from equations (19) and (20) that 

y(o) 2 = R (o) ( 21) 

and 
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( 22 ) 


R, (0) = 2 R (0) 


and integration of equation (22) gives 


R (0) = (1 - 3r e ) 


The derived functions for H (0) , R (0 \ and V (0) fulfill the boundary conditions (1 2). 
With the convenient choice of 


co(r p ) = R (0) 


the solution of the remaining equations (15) follows: 


r* 1 * = A sin r + B cos r + 2 R*°* H (1 * 


u (0 = r(o) cos r - B sin r) + 2 R (0) ^ 


t (1) = - 3 R (0)/ 2 r 2 / 2 - R (0)1 2 (A cos r - B sin t) + T (1) (r,) , 


^0) = . r(o) (. a cos r + B sin r) + 0 (1 *(r,) 


Removing secular r- terms from equations (27) and (28) gives 


</>,< °> = V (1) /R (0) . 2 R® ' H (,) 


Tj (0) = - 2 R (0)2 H (,) 


9 



Consequently, a complete base solution is dependent on H (1) , which (as well as A, B, T (1) , 
and V (,) ) may be determined from the second order solution. 

The two independent variables, using equation (24), are now 

t £ = es (31) 


t = J w(t e ) do = J R(°)' ,/2 da = -i [1 -R (0) " 2 ] . (32) 

oo e 

C. Second Order Solution 

The second order differential equations are 

wr 2 (2 > = - r/*> + u (2) 

(J^) = . ^(0 + RtO)* 1 y(2) . R (0)* 3/2 f (2) . R (0)- 2 r (l) v (l) + R (0)" S r (l) 2 

wuj® = - u, (1) + 2H (0) r (2) + 2 R (0) h (2) + 2 r (1) h (1) + R (0 > 3 ^ u (1) 

(33) 

oov 2 (2) = - V, (1) + 2 R (0) r (1) 
cuh 2 (2) = - hj (1) + v (1) 
wt 2 (2 > = - t, (1) + r (2 > 

with boundary conditions 

r (2) (0,0) = 0 

<fi (2) (0,0) = 0 (34) 

u (2) (0,0) = 0 
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v< 2) (0,0) = 0 
h (2) (0,0) = 0 
t (2) (0,0) = 0 


(34) 

(Concluded) 


To solve this set of equations we start again with the energy and v-equations 
3/2 

v (2) = 2 R (0) (- A cos r + B sin r) + V (2) (r e ) (35) 

and 

h (2) = H (2) (r e ) . (36) 

By removing improper secular r-terms from equations (35) and (36) we get a differential 
equation for H (1) 

H n (1) + 4 R (0)3 H (1) = 0 (37) 

with boundary conditions 

H (1) (r e = 0) = V (,) (r e = 0) = H^Vg-O) = 0 . 

Solving equation (37) and obeying the boundary conditions gives 

H (,) = V (1) = 0 . (38) 

This gives 

T(o) = <A (0) = 0 (39) 

Consequently, a complete zero order solution is known. Equations (1 5) will be simplified 
by removing terms which contain H (,) or V (1) . 
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A second order differential equation for the radius is derived from the first equation 
of system (33) as 

r 22 (2 > = . R(0) l/2 ri2 0) + R (0) l/2 Ua « (40) 

with 

r ]2 (,) = Aj cost - B! sin t (41) 

and 

U j (1) = R (0) ^ 2 (A! cost - Bj sin t) - R (0) (A cos t - B sin t) + 10 R (0)4 , (42) 

and the third equation from system (33) with equation (40), 

r22 ( 2 ) + r ( 2 ) = - [- 2 R (0)l/2 B, + 2 R (0)2 B] sin t - [2 R (0)l/2 A, - 2 R (0)2 A] cos r 
+ 2 R (0)2 H (2) - 8 R (0)S (43) 


The solution, from which the secular T-terms are removed, is 


r (2) = E sin t + F cos t - 8 R (0)S + 2 R (0)2 H (2) 


(44) 


The removal of the secular T-terms of the solution of equation (43) gives the following 
differential equations for the functions A (r e ) and B(t £ ): 


B, = R (0)3/2 B 
A, = R (0)3/2 A' 


(45) 


with boundary conditions 
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B(0) = 0 


( 46 ) 


A(0) = - 2 

derived from equations (25) and (26) with 
r (1) (0,0) = u (I) (0,0) = 0 . 

The solution of equations (45) and (46) gives 
A = - 2 R (0)l/2 

(47) 

B = 0 

The component of the velocity in radius direction u (2) follows as 

u (2) =1-2 R (0)2 - R (0) ^ 2 F] sin r + [R (0) E] cos r , (48) 

and 

t (2 > = [R (0)1 ^ 2 F - 4R®' 3 ] sin r - R (0 ^ 2 E cos r + ^ R (0 ^ 2 r 3 + T (2) (7 £ ) ; (49) 

the removed secular terms show that T* 1 * is dependent on H (2) as 

T, (1) = - 8 R (0)S + 2 R (0)2 H (2) (50) 

The second order term of the central angle is derivable from the second equation 
of system (33) 

V? (2) = [2 R< 0)3/2 - 4 R (0)2 - R (0) F] sin t + R (0) E cos r - R (0)3/2 sin 2r + 

(51) 


13 



and similar to equation (50) follows 


0, (1 ) = 8 R (0)1 ^ 2 + 2 R ,0) - 2 R (0)S//2 H (2) + R (0) '^ 2 V (2) . (52) 

Equations (50) and (52) show that the derivation of a complete first order solution is de- 
pendent on H (2) and V (2) . 

D. Determination of Second Order Constants 

The unknown constants of the second order solution are E, F, H (2) , V (2) , 0 (2) , and 
T^. Furthermore, H (2) and V (2 * determine <t > and T (1) from the first order solution. 

To determine these constants, it is not necessary to solve the third order equations com- 
pletely; we need only to determine the secular T-terms of this solution. 

The third order differential equations are 
co r 2 (3) = - r/ 2) + u (3) 

co ^ 2 (3) = - *>, (2) + R (0)1 v « - R(o) 3/2 r < 3 > + 2 R^ 2 r<«> r (2) - R (0) ' 2 v< 2 > r (1) 

wu 2 (3) = - u, (2) + 2H (0) r (3) + 2R (0) h (3) + 2 h (2) r (1 > + R< 0)3/2 u (2) 

+ 2 R(°W 2 r 0) U (I) (53) 

co v 2 ® — - Vi® - — R®* + 2 R® r (2) + r (1) 

2 

co h 2 (3) = - h, (2) + v« + - R< 0) ‘ l/2 u<‘> 2 

2 

co t 2 ® = - t/ 2) + r (3) 

The boundary conditions for equations (53) are 

r< 3) (0,0) = 0 
<p (3) (0,Q) = 0 
u (3) (0,0) = 0 
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(54) 



v< 3 >(0,0) = 0 
h (3) (0,0) = 0 
t (3) (0,0) = 0 


(54) 

(Concluded) 


Considering the energy equation we get with 

u d) = 4 cos 2 7 - 8 R (0 ^ ^ cos t + 4 R® 
h, (2) = H, (2) 


(55) 


the third order energy term, 


h (3) = ^ sin 2r + H (3) ( r e ) , (56) 

and collecting improper secular 7-terms yields 

H, (2) = V {2) + R (0) ' ,/2 + 2 R(° )9/2 . (57) 

Similarly, the equation for v (3) gives for the secular part, 

V/ 2) = - R (0) - 2 R (0)6 - 16 R (0)6 + 4 R (0)3 H (2) + 2 R (0) (58) 

Differentiating equation (57) gives a second order differential equation for H (2) , 

H„ (2) - 4 R (0)3 H (2) = 0 (59) 


With the boundary conditions derived from equations (34), (35), and (36), 


( 60 ) 


H (2) (0) = 0 
H, (2) (0) = - 1 
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the solution of equation (59) is 


H (2) = C R (0) 2 + K R (0) ’ l/2 


With the boundary conditions (60) we will find 



( 61 ) 


(62) 


Consequently it follows that 

H® = j [- R® > l/2 + R® 1 " 2 ] , (63) 

and from equation (54) it follows that 

V® = - j R (0) ' l/2 - i R (0)2 - 2 R(°> 9/2 - R® l/2 . (64) 

With a knowledge of H® and V® a complete first order solution is known. To determine 
E and F, derive from the radius equation of equation (53), as similarly done for r®, 

r2J <3) = . R (o) ,/2 ria ® + Uj O) R (o) l/2 (6 5) 

with 

r 12 ® = Ej cost - F, sin t , (66) 

and u 2 (3) defined by the third equation of equation (53) gives 
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r 2 2 (3) + r(3) = ' R^ 0) ^ (E'i cos t - F] sin r) 


- R (0) [- 8 R (0)?/2 + R (0) F - R (0) ' ,/2 F, ] sin t 
+ [- R (0) E + R (0) ' ,/2 E, ] cos r + 2 R (0)2 h (3) 

+ 2 R (0) h (2) r (1) + R (0)5 ^ 2 u (2) + 2 R (0)3/2 r( ° u( ° • (67) 

Collecting secular terms of r after integration of (67) gives 

Ei = R (0)3/2 E (68) 

or 

,/ 2 

E = CR (0) , (69) 

and from equation (48), using 

h (2) (0,0) = 0 , C = 0 , (70) 

we obtain 

E = 0 . (71) 

For F we will find 

F, = K {0)3 2 F + R (0)4 + 2 R (0) H (2) (72) 

with the solution using equation (63) being 

F = i R«« s/! +1-1 R®'* + R<»‘ /! !£ , (73) 

4 5 15 12 
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whereby the boundary condition, 


F(0) = 8 


( 74 ) 


is derived from equation (44) with the use of r (2) (0,0) = 0. 

From equations (50), (52), (63), and (64) we find, obeying boundary conditions, 

T (I) = y [1 - R (0)?/2 ] - ^ [1 - R (0)S ^ 2 ] + J [1 - R< 0) ] - 2 (75) 

and 

Ad) = . -i_ rCo)' 5 ^ 2 . I R (o )' 3 / 2 . I R (o) _1 . 2 R< 0 ) -1 / 2 - ]_ R ( 0 ) 3 / 2 . £ R ( 0 )' s / 2 

25 5 5 15 5 

+ R(°) 4 - in (R<°>) + L . (76) 

To derive a complete second order solution for the time t and the central angle tp, 
it would be necessary to evaluate a complete third order solution from equations (53) and 
to determine the unknown functions out of the fourth order differential equations. How- 
ever, looking at the solution so far derived, it seems reasonable to neglect the functions 
T (2) and 0 (2 * in comparison with the unbounded terms of r . T (2) and <p® are only chosen 
to fulfill the boundary conditions 


= 0 


and 


T (2) = 0 


E. Summarized Solution 

In this subsection, the derived solution is collected. The independent variables are 
calculated as 
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r e = es 


T = — — [ 1 - R (0) ] 
4e 


( 77 ) 


with 


e = — = constant 
m 


(78) 


and 


R (0) = (1 - 3 rj' 2/3 


(79) 


is the base solution for the radius of the spiral orbit. The variables of the orbit are 
approximately 


r = R (0) + e [- 2 R (0),/2 sin t] + e 2 {fcost - 8 R (0)S - j [R (0)S/2 - l]j . 
u = e [- 2 cos r + 2 R®*^ 2 ] + e 2 { [- R (0) ^ - 2 R® ] sin r | 

v = R(°)' /2 + e 2 [4 R<°> 2 cost - j R (0 >' l/2 - R (0),/2 - j R (0)2 - 2 R«» 9/2 ] 
h = - [2R (0 >]‘ 1 + e 2 ^ R (0) ' 2 - j R (0)1/2 J 

* = r + e [2 R (0) ‘ l/2 cost + 0<‘>] + e 2 { [2 K^ 2 ' 4 R( ° >2 ' R<0) F 1 sin T 
- R (0)3 ^ 2 sin 2t } 


(80) 


t = R (0)3/2 t + e [- 3 R<o )7/2 t 2 / 2 + 2 R (0) cos t + T (1) ] 


+ e 2 y- R (0 >" /2 T 3 + [R(° )l/2 F - 4 R<°> 3 ] sin t) 


where T (1) and 0 (1) are defined by equations (75) and (76) and F is defined by equation 
(73). 
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IV. AN ALTERNATE METHOD TO DETERMINE T 

In the case e = constant as considered here, we know from equation (1) that 

s 

t = J r(s) ds . (81) 

0 


Equation (81) is an independent equation. Hence, with the approximate solution for 
r[r e (s),r(s)] derived, we will be able to find t as follows: 

* 

s s 1/2 

t « J R (0) ds -2 e J R (0 ^ sin r ds . (82) 

0 0 

With equations (23), (31), and (32) we determine the derivatives 


dR (0) 

ds 


= 2e R (0) 


i/2 


(83) 


and 


dt _ 1 

Jr® ' 7Tr») s 


(84) 


and in a straightforward manner, we solve the first integral from equation (82), 

/ R (0) ds = - [1 - R (0) ~ l/2 ] . (85) 

e 

From two successive partial integrations of the second integral of equation (82), we 
determine 


s 1/2 3 

-2 e J R (0) sin t ds = 2 e [R (0 ^cosr - 1] - 4e 2 R*°* sin r + 0(e 3 ) 

0 


( 86 ) 
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or 


t = - [1 - R (0) ^ ] + 2e [R (0) sin r - 1] - 4e 2 R (0)3 sin r + 0(e 3 ) . (87) 

In numerical investigations, it was found that t determined with equation (87) is much 
more accurate than t calculated with the aid of equations (80). 

V. SOLUTION ANALYSIS 

A. Solution Characteristics 

Because of the regularizing transformation, the relationship between the velocity 
s v calculated dependent on s and the natural velocity V dependent on t is 

s v = v(s) = rv(t) = r*v (88) 


Obeying equation (88), the zero order solution derived corresponds to the standard 
circular asymptotic solution of spiral type trajectories, 


r (0) (r e ,r) • V o)2 (T e ,r) = 1 . 


(89) 


As seen from the differential equations of motion (1), the thrusting terms are of the first 
order. It follows therefrom that the zero order solution will also match the expansion solu- 
tion of the optimal thrust trajectory. The difference between tangential and optimal thrust 
appears in the first order terms which show similar characteristics; i.e., oscillations with 
slowly varying phase, frequency, and amplitude. 

The formulation of the equations of motion used causes all variables to steadily 
increase in value with superimposed oscillation. Figure 1 shows as an example for e = 0.001 
the behavior of the radius r and the energy h over the first few revolutions of the trajec- 
tory. The regularizing transformation corresponds to a scaling of the time history of the 
motion. Figure 2 shows the relationship between the regularized scaled time and the 
physically scaled time for e = 0.001. Using the energy h as an additional state variable 
simplifies the analysis of the expansion considerably. 
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ith e = 0.001 




B. Solution Limitations 


The low thrust solution derived has a singularity at 



( 90 ) 


This value corresponds to the escape condition, 


h = ’0 


and leads to 


1 

r(r € = -) -> «» 


( 91 ) 


( 92 ) 


which corresponds to the Keplerian solution for the equations of motion without thrust. 
This solution is, however, physically unacceptable; consequently, the solution presented is 
valid for energy levels less than zero only. Near escape, the thrusting acceleration and the 
gravity acceleration terms are of the same order of magnitude, and the motion ceases to 
possess two characteristic time scales. This changing in motion characteristics is because 
of the failure of the expansion presented which implies the assumption of the dominance 
of the gravity force. A complete solution of the singular escape problem would require the 
matching of the solution presented with an expansion of a different type in the vicinity of 
escape conditions. 

In addition, the expansion presented assumes that the component of the velocity in 
the radius direction u is small in comparison to the component in the circumferential direc- 
tion v. Again, near escape, this condition is not fulfilled, and the expansion solution fails 
to describe the trajectory. 


C. Numerical Comparisons 

To check the validity of the second order two-variable expansion, a comparison with 
numerically generated trajectories with different e-values was made. The terminal energy 
was always chosen as hf = 0 to explore the accuracy limits of the solution derived for each 

e-value. In general, the accuracy of the solution increases with decreasing e-values. In other 
words, for a prespecified accuracy, the validity range of the approximation will be extended 
in the direction of escape conditions with decreasing e-values. 
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Because of the high accuracy up to a certain energy level, no difference between 
both solutions would be noticeable in a diagram (such as Figures 1 and 2) showing a com- 
parison of the analytic approximation with the numerically calculated trajectory. To make 
accuracy and limitations of the two-variable expansion solutions clear, Figure 3 shows the 
validity limits for various e-values. These limits depend on the accuracy required, expressed 
in number of significance digits of the most inaccurate variable u. Figure 3 shows how 
many digits of the analytic approximation of u coincide with the numerically generated 
solution. For the numerical integration, a Runge-Kutta-Fehlberg [7] formula of the 
seventh order with stepsize control was used. 

Figure 4 shows the corresponding numbers of revolutions around the central body. 
Figures 3 and 4 show that for e < 0.001 the high accuracy of the approximate solution 
will be kept over the whole inner multirevolution part of the spiral trajectories. The loss 
of accuracy occurs in the outer part, when the variables change more rapidly. The accuracy 
in ip and t corresponds to the accuracy of u. The accuracy in r is about 1 to 2, and in 
h and v is about three to four orders of magnitude better than the accuracy of u. 

A numerically generated comparison for the example e = 0.001 is given in com- 
puter printout form in the Appendix, s is the independent variable; x = x(r, p, u, v, h, t), 
where x is the numerically generated state vector; and y = y(r, p, u, v, h, t), where y is the 
approximate solution. 


VI. CONCLUDING REMARKS 


An analytical approximate solution is presented for a low-thrust energy increase 
trajectory whereby a constant acceleration acting in a tangential direction was assumed. 

The solution was derived using a two-variable expansion method. The solution is shown to 
be very accurate for values of the disturbing acceleration e < 0.001 . This accuracy holds 
up to high energy levels, close to escape. 

The solution fails to describe the trajectory near escape because the method of 
derivation of the presented solution assumes the gravity force to be strong in comparison 
to a small disturbing force. Near escape, thrust and gravity acceleration are of the same 
order of magnitude. Furthermore, the solution implies that the velocity component in the 
circumferential direction is large compared to the component in the radius direction — a 
condition that is not valid near escape. 

Even considering those limitations, the solution presented offers a very good simpli- 
fied problem to gain experience in evaluating an accurate analytical approximation of the 
optimal low thrust trajectory. The tangential thrust trajectory solution keeps most of the 
characteristics of the solution for the optimal trajectory. The differences between those 
appear in the first order terms. 


25 



f 



THRUSTING ACCELERATION, € 


Figure 3. Energy level validity limits of the two-variable expansion solution. 
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APPENDIX 

TABULATED COMPARISON OF NUMERICALLY CALCULATED ORBITS 
AND THE TWO-VARIABLE ANALYTICAL 
APPROXIMATED SOLUTION 
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x(r, <p, u, v, h, t) NUMERICALLY GENERATED STATE VECTOR 
y(r, <p, u, v, h, t) ANALYTICAL SOLUTION 
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